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Abstract — Two rectangular wing models with a hole 
have been designed and tested in the Duke University 
wind tunnel to better understand the effects of damage. 
A rectangular hole is used to simulate damage. The wing 
with a hole is modeled structurally as a thin elastic plate 
using the finite element method. The unsteady aerody- 
namics of the plate-like wing with a hole is modeled using 
the doublet lattice method. The aeroelastic equations of 
motion are derived using Lagrange ’s equation. The flut- 
ter boundary is found using the V-g method. The hole ’s 
location effects the wing’s mass, stiffness, aerodynamics 
and therefore the aeroelastic behavior. Linear theoret- 
ical models were shown to be capable of predicting the 
critical flutter velocity and frequency as verified by wind 
tunnel tests. 

Keywords: flutter, damaged wings, finite element 
method, doublet lattice method, wind tunnel experi- 
ments 

1 Introduction 

Fighter aircraft are subject to attack by highly accu- 
rate missile defense systems when flying in enemy ter- 
ritory. Ballistic missiles can create through-hole type 
damage in the lifting surface. The aircraft may not be 
destroyed depending on the hole location and the size 
of the hole. However, the fluid-structure interaction will 
change the wing’s aeroelastic behavior. As a result, the 
flutter velocity for which the undamaged plane was de- 
signed to fly will change. 

Over the past 50 plus years the open literature is 
scarce on the aeroelastic impact of ballistic damage on 
modern aircraft. The first work in this area was per- 
formed by Biot and Arnold in 1950 [1] as reported by Dr. 
Ronald Stearman at the University of Texas- Austin [2] . 
Dr. Stearman performed extensive research on sig- 
nal detection of ballistic damage between the 1970’s 
and 1980’s. In 2007, the crack-induced effects on the 
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aeroelasticity of unswept composite wings were investi- 
gated [3]. Otherwise, the research related to ballistic 
damage in the unclassified literature has focused on in- 
vestigating the aerodynamics experimentally. In 1982, 
NASA published experimental findings on the aerody- 
namic effect of hole location caused by ballistic damage 
during certain flight conditions [4] . Other investigations 
include evaluating the effect of ballistic damage on the 
aerodynamics of helicopter rotor airfoils [5] and the ef- 
fect of the shape of the hole used to simulate damage [6] . 

Aeroelastic studies on undamaged aircraft and wind 
tunnel models are found more frequently in the liter- 
ature. Flutter and limit cycle oscillations have been 
studied experimentally and theoretically on various edi- 
tions of the actual F-16 aircraft during flight tests for 
different external store (i.e. , missiles and fuel tanks) con- 
figurations by Denegri and colleagues [7-10]. In order 
to improve the design methods for the F-16, wind tunnel 
models were designed using various theoretical models 
by Dowell, Tang, Attar [11 15] and Gordnier [16] to 
better understand the mechanisms causing flutter and 
limit cycle oscillations. 

In the present study, a rectangular wing model with 
a hole has been designed where the hole simulates dam- 
age. The wing is modeled as a thin, uniform, elastic 
plate using the finite element method. The aerodynamic 
loads are computed using the linear frequency domain 
doublet lattice aerodynamic model [17 19]. The aeroe- 
lastic model is derived by coupling the aerodynamic and 
structural models using Lagrange’s equation. 

In this paper, two rectangular models with a hole are 
designed. The hole location varied in the spanwise loca- 
tion. The theoretical flutter velocity and frequency were 
validated against tests conducted at Duke University’s 
subsonic wind tunnel. 

2 Finite Element Method 

The wing is modeled as a thin, elastic, and isotropic 
plate. These assumptions are categorized as the small 
deflection theory of bending because the “plane sections 
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initially normal to midsurface remain plane and normal 
to that surface after bending [20].” 

The finite element method uses the principle of vir- 
tual work to derive the structural model. The principle 
states that the change of strain energy due to the virtual 
displacement is equal to the change in external work due 
to applied forces, SU = 5W. The Rayleigh Ritz Method 
is used with the principle of virtual work in the commer- 
cial finite element software, ANSYS™, to formulate the 
general equations of motion for the plate-like system. 

SU = J {6e} T {a} dV (1) 

SW = J({Su} T {F} + {Su} T p{u})dV (2) 
+ j {6u} T {$}dS 

M = [JV]{d},{ii} = [JV]{d} (3) 

Equations 1, 3, and 3 represent the principle of virtual 
work on the element level used in the finite element 
method where S is the virtual operator, U is the strain 
energy, W is the external work, {«} is the displacement 
vector, {A} is the body forces vector, $ is the matrix 
of surface tractions, and p is the density. The displace- 
ment vector ({it}) is composed of nodal displacements 
(d) and shape functions (N). The nodal displacements 
are the temporally dependent nodal degrees of freedom 
and the shape functions are spatially dependent. The 
shape functions in ANSYS™ vary with the type of ele- 
ment model used. Element SHELL63, a 4-node quadri- 
lateral shell, is used here to model the uniform plate [21] . 
SHELL63 exhibits six degrees of freedom since it can 
translate and rotate in the x, y, and z directions and 
axes, respectively. The plate is clamped along the root. 
The cantilevered boundary condition is enforced by stat- 
ing that no displacement and rotation occurs where the 
plate is clamped. 

The generalized matrix form of the equations in Eqn. 
4 is formulated by ANSYS™ after substituting the 
shape functions and the assumptions into Eqns. 1-3. 
In determining the natural frequencies and modes, the 
external loads (r ext ) in Eqn. 4 are assumed to equal 
zero. The resulting equation is solved in ANSYS™ by 
further assuming simple harmonic motion to remove the 
temporal dependence. The result is the classical eigen- 
value problem that is solved using the Block Lanczos 
Method [21]. 

[' m\{d} + [k]{d} = {r (4) 

2.1 Correlation of experiment with results 
from the finite element method 


steel rectangular plate on all edges with a centralized 
hole, see Fig. 1. The analysis and experiment varies the 
hole radius, R , for a range of R/a where a=0.2032 m. 
The nondimensional frequency predicted by ANSYS™ 
is lower than those found analytically by Takahashi, see 
Fig. 1. ANSYS™ is a higher fidelity structural model 
than the one used by Takahashi. Takahaslii’s model 
used products of beam modes. 

3 Doublet Lattice Aerodynamic 
Method 

The aerodynamics for the rectangular wing with a 
hole is modeled by potential flow theory. See the aero- 
dynamic potential equation in Eqn. 5 where <f> is the ve- 
locity potential, a^o is the speed of sound, and Uoo is the 
free stream velocity. The far field boundary condition 
needed to solve Eqn. 5 states that the wave disturbances 
created by the lifting surface’s motion propagate toward 
infinity with no reflections. As a result, all the wind tun- 
nel walls are assumed to be infinitely far away accept for 
the tunnel’s floor. The floor is accounted for by using 
the method of images. The other boundary condition, 
see Eqn. 6, states the normal velocity of the fluid at 
the surface equals the normal velocity of the body. For 
the thin wing in this work, F(x,y, z,t) = z — f(x,y,t) 
in Eqn. 6 where z = 0 at the surface of the airfoil and 
f(x,y,t) refers to the height of the wing surface above 
the plane. 


V^- — 

a oo 


d 2 (j) 

~dW 


■2 U a 


d 2 (j) 

dxdt 


U; 


2 

’ dx 2 


= 0 (5) 


(JJp — » 

— + q ■ VF = 0 (6) 

The aerodynamic potential equation is transformed 
to a Kernel (integral) equation by using Green’s The- 
orem. The Kernel ( K ) equation is solved numerically 
in Eqn. 7 using a distributions of doublets to relate 
pressure (p) to downwash (w) and the doublet lattice 
method [17]. Equation 7 is rewritten in matrix form in 
Eqn. 8 where D l:] is the kernal approximated using the 
doublet lattice method. The doublet lattice method is 
a “finite element approach” to approximate the kernel. 
The doublet lattice method divides the wing into a se- 
ries of panels (boxes). One set of sides on the panels 
must be parallel to the freestream velocity. 


w(x,y ) 


1 

87T 


I< (x, £; y, ??; w, M)p(£, p)d^dp (7) 


■./ = J2 n 


j= i 


(8) 


Consider Takahashi’s experiment in the 1960’s as re- 
ported by Leissa [22]. Takahashi examined a clamped 


The doublet lattice method (DLM) is suitable to model 
the pressure across a planar wing with a hole. The DLM 
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(a) Clamped rectangular plate having a centralized circular hole 
studied by Takahashi 



(b) Results for Takahashi and ANSYS™ for a/b=l/2 and 
a/b=l 



(c) Results from Takahashi, ANSYS™ and Experiment 
a/b=l/2, v = 0.3 


methodology devised by Rodden and et al was used to 
calculate the aerodynamic influence coefficient matrix 
in a doublet lattice code written in-house [17, 19]. The 
in-house code does not use the substitution originally 
suggested by Rodden in that the steady portion cal- 
culated using the doublet lattice method should be re- 
placed with the steady portion from the vortex lattice 
method. In addition, the code uses the quartic approx- 
imation used to improve the DLM [19]. The two un- 
knowns in this problem are the pressure jump on the 
wing area around the hole and the downwash on the 
hole area. However, the latter is not needed in this 
work. In the hole region, the pressure jump equals to 
zero. The downwash on the wing portion around the 
hole is known. In determining the pressure over the por- 
tion of the wing surrounding the hole, all that is needed 
is to set the pressure to zero in the hole. The downwash 
and aerodynamic influence coefficient matrix are known, 
so its takes a simple matrix inversion to determine the 
unknown pressures. 

The local doublet strength is proportional to the pres- 
sure jump, so in the wake and the hole, the doublet 
strength is zero. The wing aerodynamics with a hole 
requires setting the doublets’ equation in the hole equal 
to zero. Computationally, it is easier to set p = 0 since 
it is the same as not having any elements there at all. 
Initially, the aerodynamic influence coefficient matrix 
(AICM) is computed for a wing without a hole using 
the doublet lattice method. The panels where the holes 
are located affects the AICM and the downwash. The 
hole on the i th panel is represented by having zeros on 
the i th row of the AICM except for on the diagonal. At 
the diagonal location a one is placed on the i th row. In 
the downwash vector a zero is placed on the i th row to 
force p = 0, see Eqn. 9. The substitutions are repeated 
for every panel on the wing that the hole region covers. 



D 12 
1 


D\n 

0 




(9) 


4 Aeroelastic Analysis 

The aeroelastic model couples the structural and the 
aerodynamic models. The aeroelastic equations are for- 
mulated using Lagrange’s equations and modal meth- 
ods. The only caveat is the structural and aerodynamic 
grids are different. 


Figure 1: Takahaslii’s experiment and a comparison us- 
ing ANSYS. Plots of nondimensional frequency pararn- 

2 4 

eter — £ p as a function of hole size for a clamped rect- 
angular plate. D = Eh 3 / 12(1 — v 2 ) is flexural rigidity. 


d dL d L 

dt dqi dqi 


(10) 


The two unknowns in Lagrange’s equation, see Eqn. 10, 
are the Lagrangian (L) and the generalized forces (Qi). 
The Lagrangian is the difference between the kinetic and 
potential energy of a structure. The generalized forces 
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are computed from the aerodynamic forces. The La- 
grangian is computed using the natural modes and fre- 
quencies obtained from the structural theory. The first 
10 out of 1000 natural modes found in ANSYS™ are 
retained for the aeroelastic model. Generally, the lowest 
modes have the greatest impact on aeroelastic behavior. 
The natural modes and frequencies are used in a modal 
series expansion. The modal series expansion, Eqn. 11, 
contains the generalized coordinate ( q m (t )) and the nat- 
ural mode (z m (x,t)). 

z a = '^2q m {t)z m (x,y) (11) 

m 

The kinetic and potential energy expressions in Eqn. 
12 and Eqn. 13, respectively requires knowing the gen- 
eralized mass (M m ) and stiffness ( K m ). The generalized 
mass is computed from Eqn. 14 where m a is the mass 
per unit area. The use of natural modes reduces the 
computation of the stiffness matrix to the mass matrix 
times the natural frequencies. 


T = 1/2^ M m q 2 m 

m 

(12) 

V=l/2 Y K mq 2 m 

m 

(13) 

M m = J J m a z^dxdy 

(14) 


The generalized aerodynamic forces are computed 
from the pressures found using the doublet lattice 
method. The downwash in Eqn. 15 is needed by the 
DLM and it comes from the natural mode deflection 
and the slope of deflection. As is, the deflection points 
of the natural mode do not match the aerodynamic grid, 
so a polynomial curve-fitting technique is applied using 
a least squares method. A fifth order polynomial is fit- 
ted to each natural mode to determine the polynomial 
coefficients needed to characterize the mode. The mode 
characterization allows the deflections to be found at 
any x and y, which is used in calculating the downwash 
for DLM. 

w _ 1 dz a (x,y,t) , dz a {x,y,t) 

u~u dt + dx ■ (15j 

The natural modes simplifies computing the mass and 
stiffness matrices for the aeroelastic model to just the 
terms along the diagonal. The polynomial’s accuracy 
is evaluated by looking at the diagonal terms in the 
mass matrix since it should approximately equal one 
when normalized. The generalized forces are computed 
from the downwash found using the polynomial in the 
doublet lattice method. Using Eqns. 11-14, Lagrange’s 
equation condenses to Eqn. 16. The aerodynamic loads 
are computed as a function of nondimensional frequency 
( k ) in the summation A mi {k)q m . 


YA m j(k)qm = Qi (16) 


The nondimensional frequency k that is needed to 
compute the aerodynamics is well suited for the V-g 
method in the aeroelastic analysis. The V-g method 
multiplies the stiffness matrix by 1 + ig where g is the 
fictitious damping adding to find the neutrally stable 
solutions. Additionally, simple harmonic motion ( q = 
(f wt ) is assumed so q can be factored out of Eqn. 16. 
Equation 16 is rearranged for the V-g method, see Eqn. 
17. At each iteration of k, oj 2 and g are calculated. 
Recalling k = tvb/U, therefore the velocity ( U ) can be 
computed since k is prescribed, to is calculated, and the 
semichord b is known. The specific velocity where flutter 
occurs corresponds to g = 0. When g equals zero the 
original aeroelastic equation (Eqn. 16) is recovered. 

j-[M] + l -^~[K\ - ^[A(k)]}{Q} = 0 (17) 

5 Results 

The two designed aeroelastic cantilever rectangular 
wings are made of polycarbonate with a thickness of 
0.001588 m, a chord of 0.1524 m (6 in), and a length 
of 0.3048 m (12 in). The plate’s material properties are 
the followimg: p=1217 kg/m 3 , E= 2.4e9 N/m 2 , and v= 
0.33. In Model A, see Fig. 2, the rectangular hole mea- 
sures with a = 0.057 m (2.25 in) and b =0.095 m (3.75 
in). The base of the hole is 0.0579 m from the root chord 
and 0.048 m from the leading edge. The hole size and 
the amount of mass removed is roughly 12% of the total 
wing. Model B is identical to Model A except the hole 
location goes outward in the spanwise location. The 
base of the hole in Model B is 0.1524 m from the root 
chord and 0.048 m from the leading edge. The struc- 



Figure 2: Dimensional drawing of rectangular wing with 
a rectangular hole near root 
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Table 1: Comparison of theory (The.) and experiment (Exp.) for first five natural frequencies (in Hz) for can- 
tilivered rectangular wing with rectangular hole near the root and tip. Bend, and Tors, refer to the bending and 
torsion mode, respectively 




No Hole 

Model A 

Model B 

mode no. 

mode char. 

The. 

Exp 

The. 

Exp. 

The. 

Exp. 

1 

Bend. 

3.99 Hz 

4.13 Hz 

3.49 Hz 

3.50 Hz 

4.20 Hz 

4.25 Hz 

2 

Tors. 

16.95 

17.25 

15.08 

15.50 

16.02 

16.13 

3 

Bend. 

24.86 

24.38 

24.69 

25.00 

23.24 

23.75 

4 

Tors. 

55.33 

54.25 

53.28 

53.75 

50.10 

51.25 

5 

Bend. /Tors. 

69.84 

69.00 

71.99 

73.65 

66.25 

68.00 


tural behavior for the rectangular plates was found the- 
oretically using the finite element method in ANSYS™. 
The structural behavior was tested experimentally us- 
ing a vibration test for the first five natural frequencies, 
see Table 1. The lower frequencies are known to have 
the largest impact on aeroelastic behavior. The struc- 
ture was excited using a Briiel & Kjaer (B&K) 4810 
mini-shaker mounted with a B&K 8200 force transducer. 
The structure’s response was recorded using the B&K 
type 4375 accelerometer placed farthest away from the 
shaker near the tip of the cantilevered plate. The shaker 
is powered by an external power amplifier that is pro- 
vided a pseudo random signal from the B&K 4 Chan- 
nel PULSE™ 3560-T-C, the front-end data acquisition 
system. PULSE™ is an advanced software/hardware 
analyzer platform developed by Britel & Kjaer for data 
acquisition and analysis. PULSE™ determines the fre- 
quency response by performing a fast fourier transform 
analysis of the acceleration time series data. 

The correlation between theory and experiment for 
the structural analysis is good. The theoretical struc- 
tural models require about 1000 finite elements to 
achieve convergence. The results are observed to be- 
have according to beam theory even though a plate is 
modeled. The natural frequencies decrease in Model A 
when the hole is closer to the root. This suggests that 
“stiffness effects” are more dominant than mass effects 
in this region. The natural frequencies increase in Model 
B when the hole is closer to the tip. This suggests that 
“mass effects” are more dominant in this region. 

The aeroelastic model was developed using the 
method outlined in the previous section. The doublet 
lattice method included 256 panels across the wing (512 
total for mirror image) with 16 divisions in the chord- 
wise direction and 16 divisions in the spanwise direction 
to calculate the generalized forces. The hole consumes 
30 panels. 

The flutter frequency was determined experimentally 
using the PULSE™ data acquisition system with an 
accelerometer fixed near the root chord on the trailing 
edge. The velocity of the windtunnel at flutter was read 
from a voltmeter connected to a pitot static tube. The 
experimental flutter velocity and frequency are generally 
in good agreement with the theory, see Table 2. 


Table 2: Theoretical and experimental flutter velocity 
and frequency listed for cantilevered rectangular wing 
with rectangular hole 



No Hole 

Model A 

Model B ! 


The. 

Exp. 

The. 

Exp. 

The. 

Exp. 

U / in m/s 

20.8 

20.05 

21.5 

20.65 

25.3 

25.2 

c Of in Hz 

10.3 

11.50 

8.5 

9.18 

8.3 

9.4 


6 Conclusions 

Flutter of a rectangular cantilever wing was predicted 
using the finite element method and the doublet lat- 
tice method in forming the aeroelastic model with La- 
grange’s equation. The flutter velocities and frequen- 
cies are reported for a rectangular wing without a hole 
and for a rectangular wing with a hole in two differ- 
ent spanwise locations. Linear models are incapable of 
determining the flutter amplitude. 

The agreement between theory and experiment for 
flutter is generally good for a hole that is 12% of the 
wing area for configurations analyzed. The interesting 
behavior is the increase in the flutter velocity from the 
case without a hole. One might anticipate the flutter 
velocity would decrease due to the presence of a hole, 
but the theory and experiment indicate the opposite. 
Further, calculations for other sizes showed expected re- 
sults. For a model with a hole size of 6% the effect of 
flutter is small. 

Developing a model for a hole size of 6% of the wing 
area yields changes that do not make a noticeable im- 
pact in investigating the stuctural behavior and the 
aerodynamics independently. Quite the opposite, a hole 
area 25% of the wing makes a large impact but the struc- 
ture is challenged from its own inertia loads before even 
considering the aerodynamics loads. 

Future work may evaluate using a time marching code 
to model the aerodynamics of a wing with a hole to de- 
velop nonlinear aeroelastic models. To date there is no 
time marching code that is capable of modeling aerody- 
namics of a wing with a hole. However, a clever substi- 
tution similar to the one implemented for the doublet 
lattice code could be potentially derived. 


RELEASED - Printed documents may be obsolete; validate prior to use. 




Acknowledgment 

The authors would like to thank Ronald Steannan of 
the University of Texas- Austin, Max Blair of Wright- 
Patterson Air Force Base, Peter Attar of the University 
of Oklahoma, and William Rodden for their insights and 
support of this work. 

References 

[1] M.A. Biot and L. Arnold, Study of vulnerability of 
aircraft to damage induced flutter, Ballistic Research 
Laboratories, 743, Oct. 1950. 

[2] G. Chen and R.O. Stearman, “A damage in- 
duced aeroelastic failure mode involving combina- 
tion and parametric resonant instabilities of lifting 
surfaces,” in A Collection of Technical Papers Part 
2: Structural Dynamics and Design Engineering, 
AIAA/ASME/ASCE/AHS 23rd Structures, Struc- 
tural Dynamics and Material Conference, May 1982. 

[3] K.H. Wang and D.J. Inman, “Crack-induced effects 
on aeroelasticity of an unswept composite wing,” 
AIAA journal, Vol. 45, No. 3, pp. 542-551, 2007. 

[4] M.L. Stearman, Wind-tunnel studies of the effects 
of simulated damage on the aerodynamics charac- 
teristics of airplanes and missiles, NASA Technical 
Manual (TM)-84588, Dec. 1982. 

[5] K.W. Robinson and J.G. Leishman, “Effects of bal- 
listic damage on the aerodynamics of helicopter ro- 
tor airfoils,” Journal of aircraft, Vol. 35, No. 5, pp. 
695-703, 1998. 

[6] P.M. Render, S. De Silva, A. Walton and M. Mani, 
“Experimental investigation into the aerodynamics 
of battle damaged airfoils,” Journal of aircraft, Vol. 
44, No. 2, pp. 539-549, 2007. 

[7] C.M. Denegri Jr., “Limit cycle oscillation flight test 
results of a fighter with external stores,” Journal of 
aircraft, Vol. 37, No. 5, pp. 761-769, Sep. -Oct. 2000. 

[8] R.W. Bunton and C.M. Denegri Jr., “Limit cycle os- 
cillation characteristics of fighter aircraft,” Journal 
of aircraft, Vol. 37, No, 5, pp. 916-918, Sep. -Oct. 
2000 . 

[9] C.M. Denegri Jr. and M.R. Johnson, “Limit cycle 
oscillation prediction using arti 

cial neural networks,” Journal of guidance, control, 
and dynamics, Vol. 24, No. 5, pp. 887-895, Sep. -Oct. 
2001 . 

[10] C.M. Denegri Jr., J.A. Dubben and D.L. Maxwell, 
“In-flight wing deformation characteristics during 
limit-cycle oscillations,” Journal of aircraft , Vol. 42, 
No. 2, p. 500-508, Mar.-Apr. 2005. 


[11] P.J. Attar, E.H. Dowell and J.R. White, “Modeling 
the LCO of a delta wing using a high fidelity struc- 
tural model,” Journal of aircraft , Vol. 45, No. 2, pp. 
1209-1217, Sep.-Oct. 2005. 

[12] D.M. Tang, J.K. Henry and E.H. Dowell, “Limit- 
cycle oscillations of delta wing models in low sub- 
sonic flow,” AIAA journal, Vol. 37, No. 11, pp. 1355- 
1362, Nov. 1999. 

[13] D.M. Tang, P.J. Attar and E.H. Dowell, “Flut- 
ter/LCO analysis and experiment for wing-store 
model, AIAA journal, Vol. 44, No. 7, pp. 1662- 1675, 
2006. 

[14] D.M. Tang and E.H. Dowell, “Flutter and limit cy- 
cle oscillations for a wing-store model with freeplay,” 
Journal of aircraft, Vol. 43, pp. 487-603, Mar. Apr. 
2006. 

[15] D.M. Tang and E.H. Dowell, “Experimen- 
tal/theoretical correlation study of gust response for 
a wing-store model with freeplay,” Journal of sound 
and vibration, Vol. 295 Nos. 3-5, pp. 659-684, 2006. 

[16] R.E. Gordnier, “Computation of limit-cycle oscil- 
lations of a delta wing,” Journal of aircraft, Vol. 40, 
pp. 1206-1208, Nov.-Dee. 2003. 

[17] E. Albano and W.P. Rodden. “A doublet-lattice 
method for calculating lift distributions on oscillat- 
ing surfaces in subsonic flows,” AIAA journal, Vol. 
7, No. 2, pp. 279-284, Feb. 1969. 

[18] W.P. Rodden, J.P. Giesing and T.P. Kalman, “Re- 
finement of the nonplanar aspects of the subsonic 
doublet-lattice lifting surface method,” Journal of 
aircraft, Vol. 9, No. 1, pp. 69-73, Jan. 1972. 

[19] W.P. Rodden, P.F. Taylor and S.C. McIntosh Jr, 
“Further refinement of the subsonic doublet-lattice 
method,” Journal of aircraft, Vol. 35, No. 5, pp. 720- 
727, Sep.-Oct. 1998. 

[20] A.C. Ugural, Stresses in plates and shell, McGraw 
Hill Companies, Inc., Boston, second edition, 1999. 

[21] Swanson Analysis Systems Inc., ANSYS user man- 
ual, release 6.1, 2002. 

[22] A.W. Leissa, Vibration of plates, NASA Technical 
Report, SP-160, 1969. 


RELEASED - Printed documents may be obsolete; validate prior to use. 



